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1 Introduction 

We proposed a novel characterization of errors for numerical weather predictions. In its 
simplest form we decompose the error into a part attributable to phase errors and a remain- 
der. The phase error is represented in the same fashion as a velocity field and is required to 
vary slowly and smoothly with position. A general distortion representation allows for the 
displacement and amplification or bias correction of forecast anomalies. 

Characterizing and decomposing forecast error in this way has two important applica- 
tions, which we term the assessment application and the objective analysis application. For 
the assessment application, our approach results in new objective measures of forecast skill 
which are more in line with subjective measures of forecast skill and which are useful in 
validating models and diagnosing their shortcomings. With regard to the objective analysis 
application, meteorological analysis schemes balance forecast error and observational error 
to obtain an optimal analysis. Presently, representations of the error covariance matrix used 
to measure the forecast error are severely limited. For the objective analysis application our 
approach will improve analyses by providing a more realistic measure of the forecast error. 
We expect, a priori , that our approach should greatly improve the utility of remotely sensed 
data which have relatively high horizontal resolution, but which are indirectly related to the 
conventional atmospheric variables. 

In this project we are initially focusing on the assessment application, restricted to a 
l'ealistic but univariate 2-dimensional situation. Specifically we study the forecast errors 
of the sea level pressure (SLP) and 500 hPa geopotential height fields for forecasts of the 
short and medium range. Since the forecasts are generated by the GEOS (Goddard Earth 
Observing System) data assimilation system with and without ERS 1 scatterometer data, 
these preliminary studies serve several purposes. They (1) provide a testbed for the use of 
the distortion representation of forecast errors, (2) act as one means of validating the GEOS 
data assimilation system and (3) help to describe the impact of the ERS 1 scatterometer 
data. 


2 Data 

Forecasts and verifying analyses made with the GEOS data assimilation and forecast system 
(Schubert et ai 1993 [11]) are used here. The particular experiments studied here are 
described by Atlas et ai (1995 [1]) in a study of the impact of ERS-1 scatterometer data 
on numerical weather prediction. The period of study is March, 1993. The forecast model 
and data assimilation system used in these experiments are identical to the GEOS-1 system 
described by Schubert et ai , except for some minor bug fixes and the modifications necessary 
to utilize surface wind vectors. Thus the control forecasts in the impact study are standard 
GEOS forecasts. In addition to the CONTROL experiment, several using different types of 
scatterometer wind information are available. Our initial prototyping and sensitivity studies 
use only the 2 x 2.5° CONTROL forecast for the period 6-11 March 1993. In addition 
we have made some comparisons to the corresponding PGLA and VARGLA forecasts. In 
all cases the CONTROL GEOS data assimilation is used as verification. The PGLA and 
VARGLA forecasts use the same setup as the CONTROL forecast, but both add ERS-1 
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scatterometer data to the CONTROL data sets in determining the initial conditions for the 
forecast. In PGLA the scatterometer data is processed using the directional filtering method 
of Offiler (1994 [8]), while in VARGLA, the variational analysis method of Hoffman (1984 
[4]) is applied to the ERS-1 cr° measurements. 


3 Methodology 

Since we require that these distortion fields vary smoothly, a spectral representation is ap- 
propriate. Determining the distortion which provides the best match is then equivalent to 
minimizing the misfit between the first field and a distortion of the second, with respect to the 
spectral coefficients of the distortion. In the present project we use a global or hemispheric 
domain, and spherical harmonics as basis functions. 

In brief, the distortion is determined by minimizing the objective function J, by varying 
the displacement and bias correction fields, where 

J = J r + Jd + Ja- 

The residual cost function, J r , measures the misfit of the distorted forecast to the verifying 
analysis. Minimizing J r improves the agreement between the (distorted) forecast and the 
analysis. The two additional penalty terms in the objective function, Jd and J a , ensure 
that the final distortion produced by the minimization is relatively smooth and not too 
large. (The terms cost function, objective function and penalty function are used more or 
less interchangeably in the literature. Here, the objective function is the quantity to be 
minimized, a cost function measures lack of fit to data and a penalty function measures lack 
of fit to a constraint.) The smoothness penalty function, Jd measures the roughness of the 
x- and y-displacements and of the bias correction, ensuring that the distortion is large scale. 
The barrier penalty function, J a , measures the magnitude of the distortion Components in a 
way so that small distortions are not penalized, but large distortions are penalized heavily. 
This has the effect of setting up a barrier to the size of the distortions which are determined. 
These last two terms are evaluated using the spectral coefficients of the distortion. 

The three terms making up J are described in the following sections. However, in our 
work so far, the spectral truncations used are so severe that Jd and J a are not used in 
obtaining the results presented here. In our studies at the beginning of the project we found 
that using the barrier and smoothness penalty functions results in distortions which are 
smaller in magnitude, but larger in scale, and residual errors which are larger in scale and 
magnitude. The spectra of the original and residual forecast error, both with and without the 
penalty functions show that a great deal of the forecast error on the scales of the distortion 
is explained by the distortion and that the penalty constraints have a strong effect limiting 
the smallest scales in the distortion. 

Note that the limits used to define J a are found to be very useful to precondition the 
minimization, even in cases where J a is not used in the functional. Simulation experiments 
demonstrated that if the control vector is scaled by its limiting values estimate, the true 
solution is quickly recovered. If the scaling is derived from the smoothing function instead, 
the minimization quickly fails with false convergence. For the case of uniform scaling of the 
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control vector, the minimization is only partially successful: the objective function is reduced 
only slowly, and after 100 iterations, only half of the original forecast error is explained. 

3.1 Residual cost function, J r 

The residual cost function J r measures the misfit between the distorted forecast and the 
verifying analysis. We denote the forecast by F, the distorted forecast by P, and the verifying 
analysis, or what is considered truth, by T. The cost function is evaluated over the global 
domain via 

UP-Tfd^ 

S.fo ' 

where the integral is the surface integral over the global domain. The distorted forecast P 
is obtained from the unmodified forecast F by adding a location-dependent bias correction 
B(\,0) to the values displaced by the displacement vector field D{\,6) = ( D U ,D V ), where 
D is expressed here in terms of its zonal and meridional components, in analogy to a wind 
field. Thus, we may write 

P(\,O) = F(\\0') + B(\,6), 

where the location (A is found by following the displacement vector D(\,0) back from 
its endpoint (A, 9). 

We represent the scalar field B by a truncated series of spherical harmonics, and the vector 
field D in terms of the spectral coefficients of the corresponding vorticity (() and divergence 
(6) fields. A degree of smoothness can thus easily be imposed by the truncation of the series, 
and further constraints can separately be imposed on the divergent and rotational parts of 
the displacement field. The control vector C for the optimization problem is thus composed 
of the spectral coefficients for £?, ( and 8: 

C = (B,(,6) t . 

* 

Both the forecast F and the verifying analysis T are available on regular latitude- 
longitude grids. For evaluation of the integral, it is convenient to first interpolate T to 
a Gaussian latitude-longitude grid, in which case the formula for J T takes the form 

^ ^ i>., - jy 2 , 

i i « 

where indices i, j denote the grid point location in longitude and latitude, Nj is the number 
of longitude points for latitude j (this number will depend on j only for reduced Gaussian 
grids), and Wj is the Gaussian weight for latitude j. These weights are normalized such that 
their sum over all latitudes is unity. 

The first step in the evaluation of the requires the spectral transformation from C 
to B tJ and (£>„, D v )ij. The next step is the evaluation of F(\',9'). Following Ritchie (1987 
[9]), we define latitude-longitude points in terms of 3-dimensional cartesian vectors centered 
on the unit sphere. The origin point (A',0'), corresponding to the 3-dimensional cartesian 
vector r, is then found in the plane of the endpoint location vector g (corresponding to 
gridpoint (A,,0j)), and the displacement vector d (corresponding to (ZA U , D v )ij): 


r = ag + ,/id, 
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where the coefficients a and (3 are chosen to satisfy the constraint that r must lie on the 
surface of the sphere, and that the length of the displacement vector d is equal to the great 
circle distance between g and r. Finally, the value F(\',0') must be obtained by bilinear 
interpolation in longitude and latitude from the surrounding grid point values. 

3.2 Smoothness penalty function, Jd 

The smoothness penalty function, Jd. is given by a simple quadratic form in terms of the 
spectral coefficients of the distortion, 

J d = E J w J C], (1) 

where j ranges over the ordering of the spectral wavenumber vectors, k, and over the compo- 
nents of the distortion — H, (, 8. Different specifications of Wj are under study. For example, 
consider the part of J d due to the bias correction. In continuous form this is given by, 

Here v is an adjustable parameter normally taken to be 1 and ctb is the scale for B. Larger 
values of v result in greater smoothing by emphasizing the contributions of higher wavenum- 
bers to J d - Using the spectral representation of B, the orthonormality of the spherical 
harmonics 'P n m , and the eigenstructure V 2 vp„ m = -n(n + we find that, 

Wk,B = {n k (n k + l)) 2l/ /<r|. 


3.3 Barrier penalty function, J a 

The functional J a serves to limit the amplitude of the distortion. For efficiency the limits 
are set on the spectral coefficients. These limits are chosen in such a way that the grid 
point (or physical space) values of bias and displacement at all locations will be limited 
by specified values. In addition the spectral coefficient limits provide a good scaling (or 
conditioning) for the minimization. Test runs using synthetic data indicate that convergence 
of the minimization is sensitive to the scaling of the control vector. 

The form of J a is chosen to be, 

J. = s ,(C,/S,) 2 ". 

where Sj are the spectral limiting values for the Cj. The adjustable parameter p, nominally 
10, controls the steepness of the barrier in spectral space (Fig. 1). 

There is no unique way of setting the spectral limits. We choose limits which correspond 
to an equipartitioning, among the spectral modes, of the contributions to the physical space 
bias correction or displacement component. Here mode means each pair (m,n). The rea- 
soning for this is that no matter what the signs of the spectral coefficients, the modes will 
tend to add up somewhere in the physical domain. On the other hand, the contributions 
within a particular mode, for example due to the sine and cosine components, are always 
out of phase and therefore add up in an rms sense. The limits on components are chosen to 
correspond to a further equipartitioning. 
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Figure 1: The component of the barrier function for a single term ( y = x 2 ^ where x = Cj/Sj), 
for ft equal to 1, 2, 10, 20. 

3.4 Implementation details 

The algorithm is implemented in Splus and Fortran. The spectral transform and computation 
of Gaussian latitudes and weights use a set of general purpose Fortran library functions. All 
computations are performed in double precision. To minimize J we use the built-in Splus 
function nlminb, which implements the algorithms of Dennis et al. (1981 [3]). The nlminb 
algorithm uses function and gradient values. Second derivatives of the cost function are 
estimated by finite differences, using repeated evaluations of the gradient and cost function. 

At first we used the version of the minimization algorithm which requires function values 
only. However, the finite difference approach is computationally inefficient, and we have 
recently developed the adjoint of the calculation of J . To develop the adjoint we use tools 
previously developed for this purpose (Hoffman et al 1992 [6]). In addition, in the present 
case the spectral transforms are nearly self-adjoint (Hoffman and Nehrkorn 1989 [7]), so 
the amount of actual adjoint code for the transforms is limited. The adjoint calculates the 
gradient of the cost function very efficiently. This technique gains more than an order of 
magnitude decrease in computational time and provides a more accurate gradient, thereby 
eliminating some difficulties reported last year. 

For example, in our first annual report, one concern was reproducibility. We reported that 
the minimization sometimes terminated immediately with a “false convergence” condition 
when starting from the reasonable initial estimate of zero distortion. As discussed, this is a 
point of maximum non-differentiability and the finite difference calculation of the gradient 
here is very sensitive to the step size. To avoid this problem, we set the initial estimate to 
random numbers. Use of random initial estimates brings with it a concern for reproducibility. 

We have made two enhancement to eliminate this problem. First, we interpolate the 
original analysis to a new grid using Gaussian latitudes and longitudes offset by half a grid 
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length. Now a zero displacement corresponds to locations interior to the grid of the analysis, 
where the interpolation of the analysis is differentiable. Second, use of the adjoint eliminates 
finite difference errors. 

Another concern was convergence. Now, with the efficiency of the adjoint code, we can 
eliminate this concern by using more than enough iterates. All results presented here use 
100 iterations. 


4 Results 

Our results show that the methodology works, that a large part of the total error may be 
explained by a distortion limited to 7T0 truncation (i.e., triangular truncation at wavenumber 
10), and that the remaining residual error contains mostly small spatial scales. With forecasts 
separated by 12 h, time continuity of the distortion fields describing the forecast errors is 
present in some regions, but lacking in other regions. 

Experiments so far are all based on the first set of forecasts described in Section 2 for 
the CONTROL, PGLA and VARGLA experiments. The CONTROL experiment analyses 
are always used as the verification. 

Initially we focused on the northern hemispheric 72 h CONTROL forecast of 500 hPa 
height valid at 00 UTC 9 March 1993. In the nominal case, we use T10 truncation and 
a Gaussian grid with half the resolution of the forecast fields to represent the distortion, 
and do not use the barrier and smoothing functions. Preliminary real data tests using T 5 
truncation were unsatisfactory and spectral analysis of the forecast and verification fields 
indicated that the minimum required wavenumber truncation is 10. Some of these results 
were reported on in our first annual technical report. More recently we have extended our 
work to the sea level pressure (SLP) field, and to the PGLA and VARGLA experiments for 
forecast ranges of 12 — 120 h. 

In general, the distortion fields appear to properly account for the forecast errors. Visual 
examination of the distortion fields reveals a number of features with phase and amplitude 
errors, some of which will be discussed below. Table 1 shows that the size of the forecast 
errors for SLP increases with time for the CONTROL forecast and that the percentage of 
the squared error explained by the distortion increases from 82% to 94%. In general, both 
the displacements and biases increase with forecast length for the first 72 h. This is obvious 
in both the plots (shown for one region below) and in statistics of the solutions. 

In Table 1 the statistics are calculated as follows for each run for the distortion corre- 
sponding to the initial estimate of zero distortion and to the solution. The rms error is the 
rms residual error for the hemisphere (in m6), calculated as the square root of J r . The rms 
gradient is the rms scaled gradient (in mb 2 ), calculated as 

where the Sj are scales given by the limiting values defined in Section 3.3. The distortion 
size is the dimensionless distance in scaled phase space, calculated as 

(B<vs,) 2 ) 1/2 - 
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Table 1: Summary statistics for SLP for the CONTROL forecast every 24 h. The calculation 
of the statistics is described in the text. 


Forecast 

Hour 

Initial rms 

Distortion 

Size 

Final rms 

Error 

Gradient 

Error 

Gradient 

24 

2.545 

3.303 

11.183 

1.094 

0.110 

48 

4.506 

7.480 

23.509 

1.787 

0.153 

72 

6.727 

10.280 

30.490 

2.199 

0.106 

96 

7.672 

12.333 

33.173 

2.295 

0.121 

120 

9.060 

17.618 

29.447 

2.217 

0.162 


In all cases 100 iterations were used by the minimization. Clearly the minimization greatly 
reduces the size of the gradient and the residual error. 


4.1 East coast cyclogenesis 

Figure 2 shows the evolution of the SLP field in the GEOS CONTROL assimilation at 24 h 
intervals beginning with the CONTROL initial conditions of 00 GMT 6 March 1993, for 
a region centered on the east coast of North America. A low forming off of New Jersey 
intensifies for 72 h and tracks ENE towards Iceland for the entire 120 h period. At 48 h a 
second low enters the domain over the Great Lakes (panel c) and follows the primary low, 
nearly merging with it at 120 h. This second low does not intensify during this period. 
Finally, a third low appears over Pennsylvania at 120 h. 

By 96 h the primary low is no longer intensifying and is now positioned somewhat ahead 
of the upper level trough at 500 mb. The 500 mb verification and CONTROL forecast are 
shown in Fig. 3, (panels a and c). The second trough, with a closed low just south of 
the Canadian maritimes is too intense and has not moved eastward enough* in the forecast. 
The main feature in the error field (forecast minus analysis) is therefore an area of negative 
difference, extending along the eastern seaboard (e). The distorted forecast (b) is much 
closer to the verification (a). The residual error (f), that is the error not explained by the 
distortion, is not only small in magnitude, but also small in scale. The distortion (d) is 
shown as vectors indicating how the forecast is displaced and as a height bias which is added 
to the displaced forecast. 

In this case, after 100 iterations, a solution for the distortion is reached which reduces 
the rms residual error from 74.29 m to 17.04 m, which corresponds to a 19-fold reduction in 
the objective function. The distortion field has significant displacements (exceeding 500 km 
in some locations) and bias corrections (100 in). 

The minimization procedure uses all degrees of freedom available to match the distorted 
forecast to the verification. Thus, the differences in the secondary trough are partly ex- 
plained by eastward displacement, and by a dipole in the bias field. In addition northward 
components of the displacement near Georgia act on the gradient of zonally averaged 500 mb 
height to weaken the trough in this area. Some of these distortion field features run counter 
to our intuitive expectations and we address possible remedies in Section 5.1.1. 

At the surface at 96 h, the secondary surface low in the CONTROL forecast (Fig. 4. a) is 
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Figure 2: Evolution of the SLP field in the GEOS CONTROL assimilation every 24 h starting 
with the initial conditions of 00 GMT 6 March 1993. 
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stronger and positioned westward relative to verification. At Boston for example, the wind 
would be stronger and with a more northerly component. The primary low is approximately 
correct in terms of strength but is postioned northwest of the verifying low by several hundred 
kilometers. The distortion (b) in this case has significant southeastward displacements in 
the area of the primary low. In the area of the secondary low eastward displacements are 
combined with positive bias corrections. The distortion for SLP is quite different from the 
distortion for 500 mb. 

The other panels in Fig. 4 show the forecasts and distortions for the PGLA and VARGLA 
experiments. In PGLA, the two forecast lows are very similar to the CONTROL case, and 
the distortions in these areas are similar. In VARGLA, the secondary low is weaker and 
more in agreement with the verification. Note that for the region plotted, the VARGLA 
distortion is larger, and the residual error is smaller. 

Interestingly, all three forecasts show evidence of a strong front in the SLP field near 
32W. In the CONTROL forecast this feature is oriented very nearly North-South, while in 
both experiments with ERS-1 data, this front curves towards the west as one moves away 
from the low. However, verification (Fig. 2.e) shows only a hint of a front in this location. 

4.2 Central Pacific Ocean 

We now focus on the region of the North Pacific Ocean near the dateline where more than one 
“best” distortion is a potential solution due to the extreme differences in features between 
the forecasts and the analyses; this results in poor time continuity of the resulting distortion 
fields. 

Figs. 5 thru 9 contain the verification, forecast, forecast error, distortion, and residual 
error fields, respectively for the 72 — 120 h (day 3-5) forecasts. In each of the figures, the left 
hand panels (a,c,e) contain the sea level pressure fields, and the right hand panels (b,d,f) 
contain those for the 500 mb height fields. 

First we compare the 500 mb height verification and forecast fields. On day 3 we see 
that the verification shows two troughs, each with closed circulations, at 40N, 153E and 
40N, 205E, respectively, separated by a broad ridge. The other important feature is another 
trough at 50-55N and 170E, just east of the Kamchatka peninsula. The corresponding 
forecast contains these features, but in a very different orientation. The easternmost trough 
is greatly underforecast in intensity, and is seen only as a weak shortwave in the flow. The 
western trough is also much weaker, and is located too far west, while the northern trough 
is too strong and is positioned too far south and east, effectively eliminating the broad ridge 
seen in the verification. By the day 4 and 5 (96 h and 120 h) time period the evolution of these 
features in the forecast and verification maps diverges greatly. The western trough intensifies 
into a major system while migrating northeastward. In response the ridge ahead of this 
system sharpens and increases in amplitude. To the east the northern and southern troughs 
coalesce at 210E. In contrast, the forecast field weakens the western trough, propagates the 
weaker eastern shortwave out of the plot region, and incorrectly amplifies and slows down 
the northern shortwave. The result is that by 120 h, the forecast and verification fields are 
essentially one-half a synoptic wavelength out of phase with one another between 165E and 
210E. This can be seen in the strong tripole pattern in the forecast error fields (Fig. 7) The 
sea level pressure fields in this case are basically a reflection of the upper level fields and 




Figure 4: Impact of ERS-1 data on the 96 h SLP field. The left column shows the CON- 
IROL, PGLA and VARGLA forecasts from top to bottom and the right column shows the 
corresponding distortions. All valid 00 GMT 10 March 1993. 
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show the same relative errors. 

The patterns of the distortion fields show only limited temporal continuity, and this may 
be due to the degree of misphasing between the forecast and verification features. This is 
particularly true for the 500 mb height field distortions (Fig. 8.b,d,f). On day 3 the distortion 
shows a large region of negative bias corrections in the center of the domain, possibly in 
response to the two large positive forecast error centers associated with the eastern and 
western troughs. With a higher spectral truncation, the algorithm may have been able to 
resolve these into two distinct features. The displacement field shows a strong northward 
component as it attempts to build the ridge seen in the verification. By days 4 and 5 the 
bias correction fields have split somewhat into two centers, but they are not sharply-defined. 
The displacement field on day 4 bears little resemblance to that on either day 3 or on day 
5. Since the features in the forecast evolve into a configuration so different from that of the 
verification there is no clearly “correct” displacement. This is consistent with our original 
assumption that the technique requires similar features in both forecast and verification 
maps in order to be effective. Nevertheless, the residual errors (Fig. 9.b,d,f) show that the 
algorithm has still reduced the forecast error. 

The sea level pressure distortions (Fig. 8.a,c,e) are somewhat more consistent in time, and 
with the features seen in the error plots. The bias corrections show two centers of negative 
values whose location and temporal evolution correspond to those in the forecast error plots. 
The displacements in days 4 and 5 are consistent with one another, with a large area of 
strong westward displacements. However, since the forecast and verification features are 
so out of phase, these displacements are not very representative of the actual displacement 
errors, even though they are again successful at reducing the forecast error, as seen in the 
residual error fields (Fig. 9.a,c,e). 

Figs. 10 and 11 show the day 5 distortions and corresponding forecast fields for the 
CONTROL, PGLA, and VARGLA experiments. In general, the distortions are qualita- 
tively similar from one experiment to another. The only exception is the sea level pressure 
displacement field for the PGLA experiment. One important characteristic of the plots in 
general is that the distortions of the pressure and height fields are often quite different from 
one another. While this is not critical for the so-called assessment application of forecast 
errors, it is important for potential use of the analysis application, when model fields may 
be changed, prior to or as part of an objective analysis procedure. 


5 Plans for future work 

We have accomplished a lot in short order, but much remains to be done. In broad terms we 
have developed and tested an efficient algorithm for determining distortions. The algorithm 
and constraints must now be tested and tuned, and applied to larger data sets. 


5.1 FY97 (the current contract) 

Our immediate plans to extend the current research fall into several related subtopics, as 
described in the following paragraphs. 
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Figure 5: Verification fields for 72 h , 96 h, 120 h (forecast days 3-5), top to bottom. Left 
column shows the verification fields for SLP and the right column shows the verification 
fields for 500 mb height. 
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Figure 6: Forecast fields for 72 h, 96 h, 120 h (forecast days 3-5), top to bottom, for the 
CONTROL experiment. Left column shows the forecast fields for SLP and the right column 
shows the forecast fields for 500 mb height. 
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Figure 7: Forecast error fields for 72 h, 96 h, 120 h (forecast days 3-5), top to bottom, for 
the CONTROL experiment. Left column shows the forecast error for SLP and the right 
column shows the forecast error for 500 mb height. 
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Figure 8: Distortion fields for 72 A, 96 A, 120 A (forecast days 3-5), top to bottom, for the 
CONTROL experiment. Left column shows the distortion for SLP and the right column 
shows the distortion for 500 mb height. 








L igure 9: Residual error fields for 72 h,, 96 h, 120 h (forecast days 3-5), top to bottom, for 
the CONTROL experiment. Left column shows residual error for SLP and the right column 
shows residual error for 500 mb height. 
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Figure 10: Distortion fields for the 120 h (day 5) forecast. Results for CONTROL, PGLA, 
and VARGLA experiments are shown from top to bottom. Left column shows results for 
SLP and right column shows results for 500 mb height. 







Figure 11: Forecast fields for the 120 h (day 5) forecast. Results for CONTROL, PGLA, 
and VARGLA experiments are shown from top to bottom. Left column shows results for 
SLP and right column shows results for 500 mb height. 
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5.1.1 More cases and variables 

The analyses will be extended to other cases in the suite of experiments. So far we have 
used only the Control and PGLA and VARGLA runs for 6 March. We have access to cases 
for 6, 11, 16 and 21 March for Control, ESA, PNMC, and PGLA experiments. A complete 
evaluation of the errors of all these experiments would constitute a major study. However 
at least the other CONTROL forecasts will be analyzed. We will also examine some of the 
Southern Hemisphere results since we know that significant impacts of scatterometer data 
exist there. 

Fields besides SLP and 500 mb height could be analyzed using distortions. Anomalies of 
SLP and 500 mb height from the zonally averaged climate mean and potential vorticity should 
provide distortions more in line with intuitive expectation. We will also experiment with 
some displacement only distortions. We will also examine differences between the current 
approach and determining first the displacement and then the bias correction, holding the 
displacement fixed. 


5.1.2 Parameter sensitivity studies 

We will examine the sensitivity of the distortion and error fields to the parameters of the 
method i/, o^, a a-, v a , p. The use of different spectral truncations and grid resolutions will be 
studied. 

Different ways of setting the spectral limits and/or scales will be explored. In particular, 
we will try to eliminate some of the assumptions in the derivation of the spectral limits by 
using the actual maximum amplitudes of the expansion functions evaluated on the transform 
grid. 


5.1.3 Refining Jd 

We will experiment with different functional forms of Jj. For the simple measures of the 
magnitude of the distortion and the lack of smoothness of the distortion, which are discussed 
by Hoffman et al. (1995 [5]), Jd takes the form of Equation 1, where the weights Wj depend 
only the total wave number n. Several different forms for Wj will be tried. For example we 
might take Wj — an 2(3 , where /? is a small integer and where a and 0 are to be determined. 
To determine the optimal values we will cross-correlate the residual errors with the distortion 
errors at various displacements. For this purpose we will vary the correlation window from 
500 to 1000 to 2000 km. If significant correlations are present we should allow greater distor- 
tions. If no significant correlations remain, we should try to further restrict the distortions. 
These ideas will be tested in parallel with synthetic data for the same versions of Jd- 

We have so far examined distortions with no constaints other than spectral truncation. 
The use of constaints like Jd and the choice of which variable to apply the distortion to 
must be explored. Differences in the distortion fields for different forecast fields should be 
examined in more detail. 
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5.1.4 Time continuity 

1 ° investigate the time continuity of the distortions we examined forecasts every 12 /?, for the 
CONTROL case of 6 March 1993. The time continuity of the distortion during the course 
of the forecast is good in some areas and involves large changes in other areas. 

Methods to enforce consistency in the evolution of the distortion must be developed. 
Consistency is most lacking if a feature in the forecast may be explained in terms of more 
than one feature in the verification. These situations are very sensitive to initial estimates 
of the distortion and the parameters of the algorithm. The use of J d and better initial 
estimates are expected to be helpful here. Also constaints like J d might be applied to 12 or 
24 h changes in the distortion. 

For example, experiments will test the suitability of extrapolating distortions in time to 
provide initial estimates for the minimization. The initial estimate of the 72 h distortion 
might be taken to be TSC^g or 2C^g f t — C 94 fo- This may also improve efficiency. 

5.1.5 Use of the 6.7 pm water vapor imagery 

The second set of our proposed experiments compare satellite data — in this case geosta- 
tionary 6.7 pm water vapor imagery — to a background field calculated from a short term 
forecast. The 6.7 pm water vapor imagery data are ideal for our study since they have strik- 
ing patterns and features, which can be matched by corresponding patterns and features 
in the short term forecast. Additionally geostationary water vapor data are available with 
high temporal frequency and near global coverage. However, we will begin our investigation 
with VAS data only, since the METEOSAT sensors do not have on-board calibration. (See 
Schmetz and Turpeinen (1988 [ 10 ]) for a discussion of the calibration of these data.) 

The approach for this task will be similar to that taken for the 500 hPa height fields. 
In this case, we take the short term forecast of the 6.7 pm water vapor imagery as Xj and 
the observed imagery as X a . There are two complications: the calculation 6 f the simulated 
6.7 pin water vapor imagery, which is discussed in the next paragraph, and the need to quality 
control the observed imagery. Quality control is required because of limited coverage, missing 
data, and the difficulty of simulating the 6.7 pm water vapor imagery at large incidence 
angles, over high terrain for dry conditions, and in the presence of cloud. It will be necessary, 
at least initially, to resample the imagery to a relatively coarse regular grid. Then we will 
determine smooth displacement and amplification fields needed to best match the forecast 
and imagery. The algorithm needed here is identical to that used for the 500 hPa height 
fields. The resulting fields of displacement and amplification provide a correction to the 
short term forecast. 

We will simulate the 6.7 pm water vapor imagery from the forecast values of temperature 
and humidity using a standard radiative tranfer model (RTM). The simulated 6.7 pm water 
vapor imagery will then be held fixed in determining the distortion. Changes in incidence 
angle related to displacements on scales of 100 km are 0(1°). The sensitivity of the calculated 
brightness temperature to incidence angle is small (Fig. 12) and will be assumed negligible 
in these calculations. The RTM used in Fig. 12, and a candidate for future calculations is 
MODTRAN. Although MODTRAN is inefficient, it is accurate and well documented (Berk 
et al. 1989 [ 2 ]). 
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Simulated GOES-7 6.7 micron Brightness Temperature 



Figure 12: Variation of 6.7 pm brightness temperature with incidence angle for the U.S. 
standard atmosphere. 
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5.2 FY98 - FYOO 

We have already submitted a proposal for a continuation of this project. We propose to 
demonstrate the usefulness of the distortion representation of forecast errors by conducting 
realistic studies of the assessment application and the objective analysis application. For the 
assessment application we will continue to study forecast errors of the 500 hPa geopotential 
height field for forecasts of the short and medium range. The forecasts and analyses which 
will be used are realistic quasi-operational forecasts. A preliminary to this task is to test and 
tune EOF representations of the distortions, using the results of the current project for the 
March 1993 forecasts. For the objective analysis application, we will apply the 6.7 pm water 
vapor imagery to determine short term forecast errors, as a standalone objective analysis 
application. These tasks are described in our proposal submitted 20 June 1997. 
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